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We model the deconfinement phase transition in quantum chromodynamics at nonzero baryon 
number density and large quark mass by extending the flux tube model (three-state, three- 
dimensional Potts model) to nonzero chemical potential. In a direct numerical simulation we conflrm 
mean-field-theory predictions that the deconfinement transition does not occur in a baryon-rich en- 
vironment. 
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0\ ■ I. INTRODUCTION 

T-H ■ 

^ ' Prospects of creating a quark-gluon plasma in the laboratory and interest in the role of such a plasma in the 

Q early Universe and in dense stars have stimulated efforts to understand the plasma/ ordinary matter phase transition 

starting from first principles in quantum chromodynamics (QCD). It is widely suspected that the phase transition 
CO , occurs at high density as well as at high temperature. While much progress has been made in characterizing the 
phase transition at high temperature at zero average baryon density in full QCD, the same degree of success has not 
CN ■ been achievable at nonzero density, chiefly because the standard SU(3) lattice action becomes complex, invalidating 
^^ ■ standard Monte Carlo methods. Simulations must then be carried out with related ensembles of real, positive weight. 
-vi ' Studies on small volumes offer intriguing hints about the phase structure at nonzero density. But with the wrong 
f— ^ . ensemble, the thermodynamic limit cannot be taken, so the phase structure cannot be ascertained [|l|,y| . Thus we turn 
f^ to simple statistical models for insight. 

T— I , The three-dimensional three-state Potts model is one of the standard paradigms for lattice QCD in the strong- 

0^ ' coupling, high-temperature, large-quark-mass limit &fl. It has been used to provide qualitative information about 
^^ \ the deconfinement phase transition. Past Potts model studies have been limited to simulations at zero quark chemical 
potential and to mean-field studies at nonzero chemical potential. In this work we show that the Potts model can be 
extended easily to nonzero baryon chemical potential, permitting direct simulation using standard techniques. 

The QCD phase transition changes character as the quark masses and flavors are varied. At zero quark mass with 
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A\ ' two or more flavors the transition restores the spontaneously broken chiral symmetry [l6[. At infinite quark mass it 

(~| , leads to "deconfinement" in the pure Yang-Mills theory. Simulations of full QCD with two flavors at zero chemical 

;iV • potential show that these two regimes are separated at intermediate quark mass by a region where no phase transition 

. ^H ! occurs — only a strong crossover. In the corresponding Potts model we find this obliteration of the deconfinement 

j^ transition, not only as the quark mass is lowered, but also as the chemical potential is increased. The former is found 

5h from direct simulation, but the latter is known until now only in mean field theory g. Here we confirm the mean 

field prediction with a direct simulation. 

Because the Potts model analogy works only in the large quark mass limit, we learn only the fate of the deconfinement 
transition and not the chiral transition. Nonetheless, to the extent that a strong crossover in QCD, observed in the 
intermediate quark-mass regime, is a vestige of deconfinement, one may speculate that a weakening of the QCD 
crossover at nonzero chemical potential is then indicated. 

In the following section we describe the model. In Sec. 3 we present results of a simulation. A summary and 
discussion is offered in the final section. 

II. THE MODEL 

We use a variant of the flux tube model of Patel H also discussed in W . The model describes a classical statistical 
system with no dynamics, consisting of a three-dimensional cubic lattice with quarks and antiquarks occupying the 
sites and color flux tubes occupying the nearest-neighbor links. A configuration is characterized by the quark number 
distribution rir € [— 3, — 2, ..., 3] and the color flux tr,i € [— IjO, 1] for each lattice site r and associated links in the 
i = 1,2,3 directions. For convenience links entering a site from a negative direction are denoted alternatively by 
Ir -7 — H-^ : „• Gauss's law is enforced in modulo 3 arithmetic: 



2.{^r,i — ^r,-i) — 'V mod 3. (1) 



i=l 



The hamiltonian assigns a mass m to each quark and flux hnk energy a to each hnk: 

if = ^a|£,,,|+^m|n,| (2) 

r,i r 

The grand canonical partition function at inverse temperature /5 and quark chemical potential fi is then 

Z/(/3,/i)= Y. exp[-P{H ^ ^^N)] (3) 

where N = '^^ Ur and the prime indicates that the sum is over aU configurations satisfying Gauss's law. 

The model is strictly static. The lowest vacuum excitations consist of mesons built from a quark-antiquark pair 
separated by one flux link and baryons consisting of three quarks on a site. Further excitations lead to extended 
baryons and mesons and more complex hadrons, always of zero triality. (We limited quark occupation to a maximum 
of three per site, but do not expect qualitative changes in our results if we increase this limit.) 

The flux tube model is equivalent to the three-state three-dimensional Potts model with complex magnetic fleld. 
That is, Zf (X Zp, where the Potts model partition function is 

Zp^J2 exP I P'J Yl R'e(2:r2;+,) + /?'/i ^ Re Zr + i(3'h' ^ Im z,. (4) 

Zr \ r^i r ^ / 

This is the form of the Potts model found from the high-temperature, high-quark-mass, strong-coupling limit of lattice 
QCD P,p[. The equivalence is established by a change of basis. The derivation starts by replacing the Gauss's law 
constraint at each site by a Kronecker delta in mod 3: 

ze2(3) 

When we introduce one such sum for each site in the lattice, the sums over link and site occupation numbers decouple 
and can be summed explicitly. The Z(3) constraint variables Zr become the "clock" spins of the Potts model. After 
rearranging the sums we have 

n (Eexp[-/3(m|nr| - M^^r)]^-"' j (6) 

r \ rij. / 

The sums under the product symbols are explicitly 

1 + 2i?e(zr<+i) exp(-/3a) (7) 

l + v!^ + v^ + z{v'^ + u) + z* [u^ + v) (8) 

where 

u = exp[— /?(m -I- /i)] and w = exp[— /3(to — /i)]. (9) 

With the aid of an identity over Z[?i) these factors can be rewritten in exponential form. The identity we need is 

1 + az + bz* = exp(c + d Re z + ie Im z) (10) 

where 



3c = ln(l +a + b)+ ln(l - a - b + a'^ + b^ - ab) (11) 

3rf = 2 ln(l + a + 6) - ln(l - a - b + a"^ + b'^ - ab) (12) 

^e^arctanff"-^)^ (13) 

Applying this identity to the first factor (0) gives the relation between the Potts spin coupling and the string energy 
(in units of the respective temperatures) 

3 V 1 - exp(— per) / 

and to the second factor (||), gives the relation between the Potts magnetic fields and the quark mass and chemical 
potential 

h(3' = I ln(l + a + 5) - ln(l - a - b + a^ + b^ - ab) (15) 
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/^'/?' = Aarctanf#i^l (16) 



where 



a = 17) 

b= ^ "^3+" 3 (18) 

Notice that h' is odd in /i and h is even. At zero chemical potential the imaginary coupling vanishes {h' = 0) and 
at infinite quark mass the real field h also vanishes. In this limit the Potts model exhibits the well known first order 
transition at a value P'^J — 0.36703(14), known from numerical simulation [g[. Here we have the Potts/flux-tube 
parameter equivalence 

(/3t J = 0.3670, P'th = P'fh' = 0) = {(3a = 1.6265, f3m = oo, /3/i = 0) (19) 

Numerical simulation has indicated that the first order transition persists to a small real magnetic field B. The 
critical endpoint is crudely known to be in the range [0.002,0.01] for which 

(/?' J = 0.365, ^'/i = [0.002,0.01],/?'/i' = 0)= {(ia = 1.632, /3m = [3.2,4.2],/?^ = 0). (20) 

Any mass m higher than this critical value should admit a first order phase transition. 

We chose to simulate the flux tube model at quark mass m/a = 5. This value was selected to assure a first order 
phase transition at zero chemical potential. Figure shows the mapping from fiux tube l3a to Potts parameters at zero 
chemical potential. We then repeated the simulation at nonzero chemical potential fi/a — 1.75. The corresponding 
map is shown in Fig. H. The Potts magnetic fields remain small. 

III. METHOD AND RESULTS 

We carry out a simulation in the occupation number basis (^, i.e. the flux tube formulation, in which the Boltzmann 
weights are real and positive and the simulation can be done easily using conventional Metropolis methods. We have 
chosen an elementary set of Metropolis moves that preserve the Gauss's law constraint, and are capable of reaching 
any valid conflguration. The local moves consist of systematically "adding" or "subtracting" one of four elementary 
color-singlet hadrons at all locations and orientations in the configuration. A single sweep of the lattice consists of 
considering each of these moves for all orientations of the hadrons at each lattice site in typewriter order. 

The Metropolis-move hadrons are these: (1) a quark and antiquark separated by one flux link, (2) a diquark and 
antidiquark separated by one fiux link, (3) a quark and diquark or antiquark and antidiquark, also separated by one 
flux link, (4) and a plaquette of flux links. The process of adding (or subtracting) a hadron consists of increasing or 
decreasing the quark occupation number and fiux link value of the configuration in mod 3 arithmetic according to 
the position and orientation of the hadron selected, respecting our exclusion principle that limits the quark number 



to the range [—3,3]. By always adding or subtracting a color singlet state, Gauss's law is always obeyed. While all 
configurations satisfying Gauss's law can be reached by a combination of these moves, this over-complete set was also 
chosen in an effort to cover the phase space efficiently. Still, we have only a local algorithm, presumably as effective 
as a local algorithm in the spin basis. For the moment we have not considered cluster algorithms analogous to those 
that have been so successful in spin systems B . 

As we have mentioned we simulate at fixed quark mass in/ a = 5 and choose two values of the chemical potential, 
namely fi/a — and 1.75. We then vary /3cr to locate the phase transition or crossover. We expect to reproduce the 
Potts first order phase transition at zero chemical potential and small real field, but at nonzero chemical potential 
we explore new territory. We simulate at a series of volumes L^ for L = 10, 20, 30, 40 in each case. In the crossover 
region we extend the simulations for, typically, 35,000 to 60,000 Metropolis sweeps. 

As expected at zero chemical potential, we find a sharp rise in energy density at f3a « 1.63 as shown in Fig. B[ To 
characterize the phase transition we study the size dependence of the peak in specific heat |TQ] . The specific heat is 
defined in the usual way in terms of the total energy E of the configuration: 

Cv/P' - HE') - {E)')/L\ (21) 

Figure H shows a peak in the specific heat at zero chemical potential that sharpens and grows with increasing 
volume. Also shown in these figures are results of a fit to a phenomenological form for a first order phase transition 



based on a simplification of the finite size analysis of Borgs et al. |11 



""'^ cosh2[7(L)(/3-/3e(L))]- ^ ' 

The simplification assumes that the probability distribution for the energy density near the phase transition receives 
a delta function contribution from the two phases with the usual Boltzmann weights: 

P{E) ex 6{E - Ed) cxp(-/3^d + Sd) + 5{E - E^) exp(-/3^„ + So) (23) 

in terms of the energy and entropy in the ordered and disordered phases. This model gives 

7 = i^(5/2 
Cy,„,ax = L'QV4 (24) 

for Q = {Ed — Eo)/L^, the latent heat of the transition. In fact the actual probability distribution resembles a sum of 
broad peaks, but since we are concerned only with extracting Cv,m&K, it suits our purpose. Although there are really 
only two independent parameters, /3c and Q, we have kept the height, width, and center unconstrained and look for 
evidence for scaling with L^. 

In Fig. |g we plot the specific heat maximum vs L^ and show that a linear relationship is plausible for L > 20. 
Errors are obtained using a bootstrap method, averaging simulation results in blocks, and extrapolating to infinite 
block size. Departures from linearity for smaller L are well known in this model and arise from finite size effects. We 
take the approximate scaling of the peak in specific heat as good evidence for a first order phase transition. 

The quark number susceptibility (three quarks = one baryon), defined as 

X, = ^-H{B')-{B)')/L' (25) 

with B = i'^rj^r/'^), the total baryon number, rises as the temperature is increased past the phase transition, as 
shown in Fig. @. 

Having tested the method, we turn to nonzero chemical potential. Shown in Fig. R is a plot of the energy density 
vs Pa. We see evidence for a weak crossover. Notice that here and in the remaining figures, we have enlarged the f3a 
scale compared with the corresponding zero chemical potential figures. Thus the crossover is far less abrupt than at 
zero chemical potential. The corresponding specific heat is shown in Fig. y. There is a diffuse peak in specific heat, 
but we see no evidence that the peak height is changing with increasing volume. The peak height is considerably 
smaller than at zero chemical potential, consistent with a broadening of the crossover. Thus we find no evidence of a 
phase transition at this chemical potential and quark mass. 

We turn to the quark number density and susceptibility. Figure ^ shows the quark number density as a function of 
f3a. We note that in the high temperature range, the quark number density at our chosen mass is quite low - namely, 
only a few per thousand sites. The density rises as temperature is increased past the crossover. As with the energy 
density, with the exception of the 10^ volume, we see no evidence for a sharpening of the crossover with increasing 
volume. The quark number susceptibility also shows no apparent trend with increasing volume, as shown in Fig. nu. 



IV. SUMMARY AND DISCUSSION 

We have extended the equivalence between the QCD-hke flux tube model and Potts model to nonzero chemical 
potential, making possible a direct numerical simulation using standard Monte Carlo methods, where only mean field 
methods have previously succeeded. We have developed numerical evidence that confirms the predictions of mean 
field theory, namely that the deconfining phase transition disappears if the baryon density is slightly nonzero. While 
our results are not applicable in the most interesting limit of small quark mass, one may speculate, nonetheless, that 
they indicate a weakening of the crossover at high baryon number density. 

A recent study by Blum, Hetrick and Toussaint analyzes the static limit of QCD and arrives at a similar conclusion 
concerning the fate of the deconfining phase transition at high density |l3]. Their resulting static action is quite 
similar to the Potts model action with imaginary magnetic field, so their method is also limited to small volumes in 
a direct Monte Carlo simulation. To obtain the comparable limit in the flux tube model, one takes m ^ oo with 
m — fj, flxed. With only two variables, namely Per and (3{m — ^), one obtains a two-dimensional subspace of the Potts 
parameter space {J(3' ,h(3\h' j3'). For large m — fj, one has /i = /i' « 0, where the flrst order deconflning transition 
occurs, whereas at small m — fj, one explores the region of large hp' and large h' (3' . In this way the results can be 
compared. 

In the corresponding canonical-ensemble study in the static limit of QCD, a weakening of the crossover signal is 
also found at nonzero baryon number |l3]. However, the authors interpret their results as suggesting the coexistence 
of two phases in the transition region. 

Could a similar change of basis also render full QCD simulations tractable at nonzero chemical potential? Normally 
the QCD path integral is formulated on a basis in which the vector potential is diagonal and the fermion integration 
is completed explicitly. The analogous change of basis diagonalizes the color electric flux and treats fermions in the 
Fock space basis, summing over states that are color singlets in quark and gluon content. At finite Nc it appears that 
the complexities of enforcing color neutrality in such a basis are prohibitive. 
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FIG. 1. Mapping of flux tube parameters f3cr to real Potts magnetic fleld P'h at fixed quark mass m/a — 5 and zero chemical 
potential. A first order phase transition is expected at /3cr « 1.63. 
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FIG. 2. Mapping of flux tube parameters /3a to Potts real and imaginary magnetic fields f3'h and f3'h' at fixed quark mass 
m/a = 5 and nonzero chemical potential fi/a — 1.75. 
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FIG. 3. Energy density vs 13a in the critical region at zero chemical potential for various lattice sizes: L — 10, crosses; 20, 
diamonds; 30, octagons; 40, squares. 
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FIG. 4. Specific iieat vs /3a at zero chemical potential for various lattice sizes. Symbols as in Fig. 
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FIG. 5. Peak in specific heat vs lattice volume at zero chemical potential 
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FIG. 6. Quark number susceptibility vs /3cr at zero chemical potential for various lattice sizes. Symbols as in Fig. 
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FIG. 7. Energy density vs Pa at chemical potential /3/i = 1.75 for various lattice sizes. Symbols as in Fig. ra. 
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FIG. 8. Specific heat vs /3a at chemical potential /3/x = 1.75 for various lattice sizes. Symbols as in Fig. H. 
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FIG. 9. Quark number density vs jSa at chemical potential /3/i = 1.75 for various lattice sizes. Symbols as in Fig. 
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FIG. 10. Quark number susceptibility vs j3a at chemical potential /3/i = 1.75 for various lattice sizes. Symbols as in Fig. H. 
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